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ABSTRACT 

Plasma instabilities (e.g., Buneman, Weibel and other two-stream instabilities) excited in collision- 
less shocks are responsible for particle (electron, positron, and ion) acceleration. Using a new 3-D 
relativistic particle-in-cell code, we have investigated the particle acceleration and shock structure 
associated with an unmagnetized relativistic electron-positron jet propagating into an unmagnetized 
electron-positron plasma. The simulation has been performed using a long simulation system in order 
to study the nonlinear stages of the Weibel instability, the particle acceleration mechanism, and the 
shock structure. Cold jet electrons are thermalized and slowed while the ambient electrons are swept 
up to create a partially developed hydrodynamic (HD) like shock structure. In the leading shock, 
electron density increases by a factor of < 3.5 in the simulation frame. Strong electromagnetic fields 
are generated in the trailing shock and provide an emission site. We discuss the possible implication 
of our simulation results within the AGN and GRB context. 

Subject headings: relativistic jets: Weibel instability - shock formation - electron-positron plasma, 
particle acceleration, magnetic field generation - particle-in-cell 



1. INTRODUCTION 

Particle-in-cell (PIC) simulations can shed light on the 
microphysics within relativistic shocks. Recent PIC sim- 
ulations show that particle acceleration occurs within the 
downstream jet (e.g., Frederiksen et al. 2004; Nishikawa 
et al. 2003, 2005, 2006, 2008, 2009; Hededal et al. 2004; 
Hededal & Nishikawa 2005; Silva et al. 2003; Jaroschek 
et al. 2005; Chang, Spitkovsky & Arons 2008; Dieck- 
mann, Shukla, & Drury 2008; Spitkovsky 2008a,b; Mar- 
tins et al. 2009). In general, these simulations confirm 
that a relativistic shock in weakly or non magnetized 
plasma is dominated by the Weibel instability (Weibel 
1959). The associated current filaments and magnetic 
fields (e.g., Medvedev & Loeb 1999) accelerate electrons 
(e.g., Nishikawa et al. 2006) and cosmic rays, which 
affect the pre-shock medium (Medvedev & Zakutnyaya 
2009). 

In this paper we present new three-dimensional simu- 
lation results for an electron-positron jet injected into an 
electron-positron plasma using a long simulation grid. A 
leading and trailing shock system develops with strong 
electromagnetic fields accompanying the trailing shock. 

2. SIMULATION SETUP 
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The code used in this study is an MPI-based par- 
allel version of the relativistic electromagnetic particle 
(REMP) code TRISTAN (Buneman 1993; Nishikawa 
et al. 2003, Niemiec et al. 2008). The simulations 
have been performed using a grid with [L^, Ly, L^) — 
(4005,131,131) cells and a total of 1 billion parti- 
cles (12 particles/cell/species for the ambient plasma) 
in the active grid. The electron skin depth. As — 
c/ujpc = 10. OA, where ujpe = (e^na/eoWe)^^^ is the elec- 
tron plasma frequency and the electron Debye length 
Ad is half of the cell size, A. This computational do- 
main is six times longer than in our previous simula- 
tions (Nishikawa et al. 2006; Ramirez-Ruiz, Nishikawa & 
Hededal 2007). The jet-electron number density in the 
simulation reference frame is 0.676 rin, where na is the 
ambient electron density, and the jet Lorentz factor is 
= 15. The jet-electron/positron thermal velocity is 
^^j,th = 0.014 c in the jet reference frame, where c = 1 
is the speed of light. The electron/positron thermal ve- 
locity in the ambient plasma is t^a.th = 0.05 c. As in our 
previous work (e.g., Nishikawa et al. 2006) the jet is in- 
jected in a plane across the computational grid located at 
X — 25A in order to eliminate effects associated with the 
boundary at a; = Xmin- Radiating boundary conditions 



are used on the planes at x 



and X 



and pe- 



riodic boundary conditions on all transverse boundaries 
(Buneman 1993). 

The jet makes contact with the ambient plasma at a 
2D interface spanning the computational domain. Here 
the formation and dynamics of a small portion of a much 
larger shock are studied in a spatial and temporal way 
that includes the spatial development of nonlinear satu- 
ration and dissipation from the injection point to the jet 
front defined by the fastest moving jet particles. 

3. SIMULATION RESULTS 

Figure la & b show the averaged (in the y — z plane) 
(a) jet (red), ambient (blue), and total (black) electron 
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density and (b) electromagnetic field energy divided by 
the total jet kinetic energy {El = X]i=e p "^i*^^ (7j ~ 1)) 
at t = 3250 Wpe^. Here, "e" and "p" denote electron and 
positron. Positron density profiles are similar to electron 
profiles. Ambient particles become swept up after jet 
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Fig. 1. — Averaged values of (a): jot (red), ambient (blue), 
and total (black) electron density, and (b): electric (red) and 
magnetic (blue) field energy divided by the jet kinetic energy at 
t = 3250 uipB ■ Panel (c) shows the evolution of the total electron 
density in time intervals of St = 250 Wpe ■ Diagonal lines indicate 
motion of the jet front (blue; < c), predicted contact discontinuity 
speed (green: ~ 0.76 c), and trailing density jump (red: ~ 0.56 c). 

electrons pass x/A ^ 500. By t = 3250 uj~^, the density 
has evolved into a two-step plateau behind the jet front. 
The maximum density in this shocked region is about 
three times the initial ambient density. The jet-particle 
density remains nearly constant up to near the jet front. 

Current filaments and strong electromagnetic fields ac- 
company growth of the Weibel instability in the trail- 
ing shock region. The electromagnetic fields are about 
four times larger than that seen previously using a much 
shorter grid system {L^ = 640 A). At t = 3250 uj^^, 
the electromagnetic fields are largest at x/A 1700, 



and decline by about one order of magnitude beyond 
x/A — 2300 in the shocked region (Nishikawa 2006; 
Ramirez-Ruiz, Nishikawa & Hededal 2007). 

Figure Ic shows the total electron density plotted at 
time intervals of St = 250 Wp^f. The jet front propa- 
gates with the initial jet speed (< c). Sharp RMHD- 
simulation shock surfaces are not created (e.g., Mizuno 
et al. 2009). A leading shock region (linear density in- 
crease) moves with a speed between the fastest moving 
jet particles < c and a predicted contact discontinuity 
speed of ~ 0.76 c (see §4). A contact-discontinuity re- 
gion consisting of mixed ambient and jet particles moves 
at a speed between ~ 0.76 c and the trailing density jump 
speed 0.56 c. A trailing shock region moves with speed 
< 0.56 c, note the modest density increase just behind 
the large trailing density jump. 
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Fig. 2. — Phase-space distribution of jet (red) and ambient (blue) 
electrons at t = 3250 LUp^ ■ About 18,600 electrons of both species 
are selected randomly. 

Figure 2 shows the phase-space distribution of jet (red) 
and ambient (blue) electrons at t = 3250 uj^^ and con- 
firms our shock-structure interpretation. The electrons 
injected with jjVx ~ 15 become thermalized due to 
Weibel instabililty-induced interactions. The swept-up 
ambient electrons (blue) are heated by interaction with 
jet electrons. Some ambient electrons are strongly accel- 
erated. 

^'lO 




Fig. 3. — Velocity distributions at t = 3250 ujpc . All jet (red) and 
all ambient (blue), and at x/A > 2300 jet (orange) and ambient 
(green) electrons are also plotted. The small (red) peak indicates 
jet electrons injected at 7^ = 15. 

Figure 3 shows the velocity distribution of all jet and 
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ambient electrons in the simulation frame. The small 
peak indicates electrons injected at 7j = 15. Jet electrons 
are accelerated to a non-thermal distribution. Ambient 
electrons are also accelerated to speeds above the jet in- 
jection velocity. The velocity distributions of jet and 
ambient electrons near the jet front (at a;/ A > 2300) are 
also plotted. The fastest jet electrons, 7 > 20, arc lo- 
cated near the jet front. On the other hand, the fastest 
ambient electrons are located farther behind the jet front 
(at a;/A < 2300). Thus, strong acceleration of the am- 
bient electrons accompanies the strong fields associated 
with the Weibel instability. 

4. DISCUSSION 

Our collisionless-shock structure can be compared to 
1-D hydrodynamic (HD) shock predictions (e.g.. Bland- 
ford & McKee 1976; Zhang & Kobayashi 2005). The 
speed of the contact discontinuity (CD) is given by ram 
pressure balance in the CD frame. Our initial con- 
ditions allow us to set the total energy density e = 
pc^ -|-p/(r — 1) = pc^ and pressure p = 0, so that the 
speed in the ambient frame becomes (Rosen et al. 1999) 



/?cd = [(7j'?^/')/(7j'?'/' + l)]/?j, 



(1) 



where r] = P]/ Pa.{= "nicn-Jmcng) and mass densities are 
determined in the "jet" and "ambient" proper frames. 
In the simulation nj = 0.0451na and 7j = 15, and fi^d = 
0.759 (7cd = 1-54) is the predicted CD speed. Formally 
this should represent the average speed of particles in the 
CD region. 

The leading shock moves at a speed given by 



2 ^ (7cd + l)[rsa(7cd -!) + !]' 

r3a(2-r,a)(7cd-i) + 2 



(2) 



where 5/3 > Fga > 4/3 is the shocked ambient adiabatic 
index. Thus the leading shock speed is predicted to be 
0.865 > > 0.783 (2 > 71^ > 1.6) where upper and 
lower limits correspond to upper and lower limits of Fga, 
respectively. 
The jump condition at the leading shock is 



(3) 



where riga is the shocked ambient density in the proper 
(CD) frame and we find 5.34 < Usa < 9.15 ria, where 
the lower and upper limits correspond to the upper and 
lower limits to Fga, respectively. Measured in the ambient 
(simulation) frame the shocked ambient density should 
be 8.2 ria < 7cd'^sa < 14.1 ria- Formally this should 
represent the total density of particles in the shocked- 
ambient region. 

Computations associated with the trailing shock are 
most easily performed in the jet rest frame designated 
below as the "primed" frame. In this frame the CD 

moves with speed = -/3cd)/(l-/?j/3cd) = -0.984 
and 7^ J = 5.60. The speed of the trailing shock in the 
jet frame, 7(3 is given by eq. (2) but with 7cd —>■ 7cd ^^"^ 
Tsa ^ Tsj where Fsj is the shocked-jet adiabatic index. In 
the jet frame 10.4 > 7^^ > 7.4 and 0.995 > -(3[^ > 0.991, 
where upper and lower limits correspond to upper Tgj = 
5/3 and lower Fsj =4/3 limits to Fgj, respectively. The 
trailing shock speed in the ambient (simulation) frame 
is 0.35 < As = (/3j - AJ/(1 - < 0.61 where the 



lower and upper limits correspond to the upper and lower 
limits of Fsj , respectively. 



The density jump at the trailing shock is given by 
eq. (3) but with 7cd — > 7cd ^^'^ ^sa — * ^sj where 
now risa/na — >■ nsj/rij where nj = 0.0451 ria with re- 
sult that the proper density of shocked jet material is 
0.70 Ha < Ugj < 1.15 ria where lower and upper limits 
correspond to upper and lower limits to Fgj , respectively. 
In the ambient (simulation) frame the shocked jet den- 
sity should be 1.08 ria < 7cd?^sj < 1'76 ria- Formally 
this should represent the total density of particles in the 
shocked jet region. 

In the simulation the speed of the trailing density jump 
is ~ 0.56 c, which is in the predicted range 0.35 < /?ts < 
0.61, a typical speed within the density-plateau region, ~ 
0.75 c, is close to /3cd = 0.76. The poorly defined leading 
shock structure moves at a speed between ^ 0.76 c and 
< c, consistent with the predicted 0.78 < P\s < 0.86. 

In the simulation the maximum density increase ob- 
served in the ambient (simulation) frame is 7cd?isa/'"'a ~ 
3.5 behind the leading shock (see Fig. la). This is 
about a factor of ~ 3 smaller than the predicted increase, 
8.2 < 7cdnsa/na < 14.1, for a fully-developed leading 
shock. On the other hand, the density increase observed 
in the ambient (simulation) frame of 7cdnsj/na }t 1 just 
before the trailing large density jump is comparable to 
that predicted, 1.08 < 7cd'^sj/''^a < 1.76, for a fully de- 
veloped trailing shock. 

Our present results can be compared to those found 
in the 2-D simulations of Chang et al. (2008) (see also 
Spitkovsky 2008a). Their simulations were performed 
in the CD frame, and material with proper density, n, 
moved into the contact discontinuity with a Lorentz fac- 
tor 7 = 15. A shock moved away from the CD with the 
predicted speed 



/?s = (Fs - 1) 



7- 



1 



-,1/2 



7 + 1 



0.47 , 



and predicted density jump 

Th ^ 1 rs7 -hi 
7n 7 Fg - 1 



3.13 



(4) 



(5) 



for a shocked adiabatic index of Fg = 3/2. 

In our simulation we have two shocks that move away 
from the CD. For our leading shock, the ambient plasma 
moves relative to the CD at a speed equal to /3cd = 0.759 
and 7 = 7cd = 1-54 in eqs. 4 & 5. In the CD frame (3s = 
0.23 and the observed density jump becomes nsa/7cdna = 
4.3 for Fg = 3/2. So we see that our leading shock speed 
would be about 50% less than that in Chang et al. (2008) 
and our density increase would be about 50% larger for a 
fully-developed leading shock in the CD frame. For the 
trailing shock, the jet moves toward the CD at a speed 
equal to -(3'^^ = 0.984 and 7 = 7^^ = 5.60 in eqs. 4 k 
5. In the CD frame f3s = 0.417 and the observed density 



increase becomes 



si/7cd"j 



3.36 for Fg = 3/2. So we 



see that our trailing shock speed would be about 11% less 
than that in Chang et al. (2008) and our density increase 
would be about 7% larger for the fully developed trailing 
shock in the CD frame. The parameters associated with 
our trailing shock are similar to those found in Chang et 
al. (2008), and the Weibel filamentation structures are 
comparable but now studied in full 3-D. 
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5. CONCLUSION 

The present simulation finds for the first time a rel- 
ativistic shock system comparable to a predicted rela- 
tivistic HD shock system consisting of leading and trail- 
ing shocks separated by a contact discontinuity, albeit 
not yet fully developed. One remarkable aspect of this 
shock system lies in the generation of large electromag- 
netic fields, up to 30% of the kinetic energy density, as- 
sociated with the trailing shock. Electromagnetic fields 
in the leading shock and contact-discontinuity region are 
over one order of magnitude lower. The large value for 
SB 0.3 in our trailing shock hints that Poynting-flux- 
dominated ejecta may not be required to explain some 
GRB observations (McMahon et al. 2006). 

Visualization of our dual shock system in the ambi- 
ent (simulation) frame provides a picture of the shock 
structure that should exist at the head of a relativis- 
tic astrophysical jet, 7jt = 15, that is less dense than 
the surrounding medium, rijt/fiam = 0.045. Within the 
AGN context, here we identify our trailing shock with 
the "jet" shock that decelerates the relativistic jet and we 
would expect synchrotron emission to originate from the 
strongly magnetized structure. Little synchrotron emis- 
sion would originate from the weakly magnetized "bow" 
shock in front of the contact discontinuity. This in fact 
is what is observed at the leading edge of extra-galactic 
jets where synchrotron emission from the bow shock is 
not typically observed. 

Visualization of our dual shock system in the "jet" 
frame provides a picture of the shock structure that 
would accompany a relativistic blast wave driven by 
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emission from this strongly magnetized structure. Lit- 
tle synchrotron emission would originate from the low 
Lorentz factor, weakly-magnetized "reverse" shock mov- 
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Our present simulation involves an electron-positron 
jet and ambient medium. We might expect similar shock- 
structure development in electron-ion simulations, albeit 
on much longer temporal and spatial scales. 
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